from mpmath import *
import math
import random




def system_matrix(delta0, epsi, delta0_deviation, epsi_deviation):
	while 1:
		delta_mod = random.gauss(delta0, delta0_deviation)
		epsi_mod = random.gauss(epsi, epsi_deviation)
		if delta0 > 0 and epsi > 0:
			return matrix([[0, 1], [1, 0]]) * delta_mod + matrix([[1, 0], [0, -1]]) * epsi_mod

def system_extra_matrix(delta0_extra, epsi_extra, t = 0, stepsize = 0, t0 = -1):
	if math.fabs(t - t0) < 0.5 * stepsize and t0 > 0:
		return matrix([[0, 1], [1, 0]]) * delta0_extra + matrix([[1, 0], [0, -1]]) * epsi_extra
	else:
		return 0
	


def system_initial_state(init):
	if init == 0:
		return matrix([1, 0, 0, 0])
	if init == 1:
		return matrix([0.5, 0.5, 0.5, 0.5])

	
